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Abstract 

We study the bifurcations of a set of nine nonhnear ordinary differential equations that de- 
scribe the regulation of the cyclin-dependent kinase that triggers DNA synthesis and mitosis in the 
budding yeast, Saccharomyces cerevisiae. We show that Clb2-dependent kinase exhibits bistability 
(stable steady states of high or low kinase activity). The transition from low to high Clb2-dependent 
kinase activity is driven by transient activation of Cln2-dependent kinase, and the reverse transition 
is driven by transient activation of the Clb2 degradation machinery. We show that a four- variable 
model retains the main features of the nine-variable model. In a three-variable model exhibiting 
birhythmicity (two stable oscillatory states), we explore possible effects of extrinsic fluctuations on 
cell cycle progression. 
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I. INTRODUCTION 



The cell cycle is the sequence of events by which a growing cell replicates all its compo- 
nents and divides them evenly between two daughter cells . Many theoreticians have 
understood the cell cycle as a periodic process driven by a biochemical limit cycle oscillator 
. However, a growing body of experimental and theoretical evidences indicates that 
the eukaryotic cell cycle is a toggle switch between two stable steady states, controlled by 
checkpoints 

mm 

\v\ . This point of view was adopted by Chen et al. in a recent 
mathematical model of the budding yeast cell cycle, and bistability in the yeast ce ll cy cle 
control system has been confirmed recently by experiments from Cross's laboratory 10|. 

Bifurcation theory is a mathematical tool for characterizing st eady state and oscillatory 
solutions of a system of nonlinear differential equations (ODE's) |l2l. Il3l|. The goal of this 
work is a detailed bifurcation analysis of Chen's model. Our bifurcation analysis supple- 
ments the numerical simulation carried out b;^_Chen et al. and clarifies their quantitative 
comparisons between experiment and theory [10]. In addition, bifurcation theory helps us 
to identify control modules within Chen's complicated model, thereby bringing some new 
insights to the yeast cell cycle control mechanism. A more through understanding of cell 
cycle control in yeast can be very helpful in future efforts to model mammalian cell cycle 
controls 

This paper is organized as follows. In Section II, we give a brief introduction to the 
budding yeast cell cycle. In Section III, we introduce Chen's model and present its one- 
parameter bifurcation diagram. In Section IV, we study saddle node bifurcations in Chen's 
model in order to provide a rigorous foundation for interpreting Cross's experiment on 
bistability of the control system jl^. In Section V, we propose a reduced model with four 
time-dependent variables, which retains the main dynamical characteristics of the extended 
model. We characterize this model using two-parameter bifurcation diagrams. In Section 
VI, we further reduce Chen's model to three variables and demonstrate that the abbreviated 
.odel d.play= to ™ complex .ode. g H Q. 

In Section VII, we use the three-variable model to study effects of extrinsic fluctuations. 
The closing section is devoted to discussion. The nine-variable mathematical model and its 
parameters are given in the Appendix. 



2 



II. A BRIEF INTRODUCTION TO THE BUDDING YEAST CELL CYCLE 



Cell cycle phases. The cell cycle is the process by which one cell becomes two. The most 
important events in the cell cycle are replication of the cell's DNA and separation of the 
replicated DNA molecules to the daughter cells. In eukaryotic cells, these events (replication 
and separation) occur in temporarily distinct stages (S phase and M phase, respectively). S 
and M phase are separated in time by gaps called Gl and G2 phases. 

During S phase ("synthesis"), double-stranded DNA molecules are replicated to produce 
pairs of sister chromatids. During M phase ("mitosis"), sister chromatids are separated so 
that each daughter cell receives a copy of each chromosome. The Gl checkpoint mechanism 
controls the initiation of S phase, and a G2 checkpoint mechanism controls entry in M 
phase. A mitotic checkpoint controls the transition from M phase back to Gl phase. The 
checkpoints monitor cell size, DNA damage and repair, DNA replication, and chromosome 
alignment on the mitotic spindle. 

Molecular controls of budding yeast cell cycle. Based on current knowledge about the 
molecular components controlling progression through the budding yeast cell cycle, a molec- 
ular wiring diagram was proposed by Chen et al. A slightly simplified version of their 
diagram is presented in Figure 1. The molecular components can be divided into four groups: 
cyclins, inhibitors, transcription factors, and proteolytic machinery. 



There are two families of cyclins in Figure 1: Gin's and Gib's p6^. These cyclins combine 
with kinase subunits (Gdc28) to form active cyclin-dependent kinase heterodimers that 
trigger cell cycle events (Gdc28/Gln2 initiates budding, Gdc28/Glb5 initiates DNA synthesis, 
Cdc28/Glb2 initiates mitosis). Gdc28 subunits are in constant, high abundance throughout 
cell cycle; hence, the activity of Gdc28/cyclin heterodimers is controlled by the availability 
of cyclin subunits. For this reason, Gdc28 is not shown in Figure 1; only the cyclin subunits 
are specified. (Each cyclin molecule is understood to have a Gdc28 partner.) 

Sicl (in Figure 1) is a cyclin-dependent kinase inhibitor: it binds to Gdc28/Glb dimers to 
form inactive trimers (Gdc28/Glb/Sicl). Sicl does not bind to or inhibit Gdc28/Gln dimers. 

Mcml, MBF, SBF and Swi5 are transcription factors for synthesis of Glb2, Glb5, Gln2 
and Sicl, respectively. 

The degradation of these proteins is regulated by a ubiquitination pathway. Proteins 
destined for degradation are first labeled by attachment of multiple ubiquitin molecules. 
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Ubiquitin moieties are attached to Clb2 and Clb5 by the APC (anaphase promoting complex) 
in conjunction with either Cdc20 or Hctl. Sicl is ubiquitinated by a different mechanism(the 
"SCF"), which (unhke the APC) requires that its substrates be phosphorylated. 

Budding yeast cells progress through the division cycle as the levels of the species in 
Figure 1 come and go. Thus the problem of cell cycle control is to understand the temporal 
fluctuations of these species. Because the species in Figure 1 are directly or indirectly inter- 
acting with all other species, simultaneous determination of their fluctuating concentrations 
require a precise mathematical model. Using mass action and Michaelis-Menten rate laws, 
the complex wiring diagram in Figure 1 can be converted into ordinary differential equations, 
and from them the molecular levels can be computed Q]. 



III. A BIFURCATION DIAGRAM OF CHEN'S MODEL 

The model proposed by Chen et al. ^ includes about a dozen ODE's and eleven algebraic 
equations with more than 50 parameters. (Refer to ^ for a complete description of the 
wiring diagram and a derivation of the mathematical model, as well as for estimates of the 
rate constants in the model.) 

In the appendix we present a reduced version of Chen's model to be used in this paper 
for bifurcation analysis. From the original model, we drop the target variables (spindle and 
bud formation, and DNA synthesis) because they are decoupled from the rest of the model. 
We reserve mass as the principal bifurcation parameter. We use the same parameter values 
as Ref . P , and they are presented in t^ appendix. Table I. 

Using the software package AUTO 17], we created a one-parameter bifurcation diagram 
(Figure 2) of the budding yeast cell cycle model, Eqn. (A1-A20), for parameter values 
given in Table I. Two saddle node bifurcations, at M 0.97 and M ^ 0.6, connect the 
stable steady states in Figure 2. There is also a subcritical Hopf bifurcation on the upper 
branch of steady states at M ^ 0.82 from which a branch of unstable limit cycles originates. 
These unstable oscillations disappear at an infinite-period saddle loop (SL) bifurcation near 
M ^ 0.73. A second branch of limit cycle oscillations, shown by filled circles that disappear 
at a different SL bifurcation point (M ^ 0.78), are stable. 

The stable steady states (solid line) at values of [Clb2]< 10"'^ represent the Gl phase of 
the cell cycle. The stable oscillatory states (filled circles) represent autonomous progression 
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through S, G2 and M phases of the cell cycle and then back into S phase. 

To get a full picture of cell cycle events, we must combine the dynamics of the cyclin- 
dependent kinase "engine" (as summarized in Figure 2) with equations for cell growth and 
division (changes in cell mass, M). To this end, we supplemented Eqn. (A1-A20) with 
an equation for mass growth, M{t) = M{0)e^^, or, in differential form, M = fiM, and a 
rule for cell division (M reset to JM whenever Clb2-dependent kinase activity drops below 
0.005). Following Chen et al. '0\ we choose / = 0.0043 because budding yeast cells divide 
asymmetrically. 

With these changes, we compute a solution of the full system, Eqn. (A1-A20) plus the 
dynamics of M, and plot the resulting "trajectory of motion" on the bifurcation diagram 
(the red line in Figure 2). This trajectory shows that the control system stays in the Gl 
phase if M < 0.97. As M increases further, the control system is captured by the stable 
limit cycle. As a result, [Clb2]T increases abruptly, driving the cell through S phase into M 
phase, then [Clb2]T drops below 0.005, causing the cell to divide and the control system to 
return to the stable Gl state. 

This bifurcation diagram of the Chen et al. model exhibits the same features of cell cycle 
models of frog eggs ,15] and fission yeast 16], namely, saddle node bifurcations associated 
with stable and unstable oscillations. Yet, there are subtle differences in these bifurcation 
diagrams. In the frog egg and fission yeast models, the large amplitude stable limit cycles 
end at a saddle-node invariant-circle (SNIC) bifurcation, not a SL bifurcation. In our case 
stable oscillations coexist with the stable steady states over a small range of mass values 
(0.78 < M < 0.95). However, when the budding yeast cell cycle model is supplemented by 
the mass growth equation, such differences seem unimportant. 



IV. SADDLE-NODE BIFURCATIONS DRIVEN BY CLN2 AND CDC20 

Recently, Cross et al. Q] experimentally confirmed bistablity in activity of Clb2- 
dependent kinase in budding yeast cells. It is interesting to mention that this result was 
predicted by a schematic sketch (Figure 9 of Ref. 0]) intuitively drawn from interrelations 
of Cdc28/Clb2 with the Gl phase cyclin Cln2 and the APC specificity factor Cdc20. We 
confirm this informal prediction of Chen et al. by a rigorous bifurcation analysis of their 
model. 
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In their experimental work, Cross et al. constructed a strain that under different experi- 
mental conditions may lack activities of Cln2 or Cdc20 or both. (To be precise, Cross et al. 
used Cln3 in place of Cln2, but that technical detail makes no difference to our analysis.) By 
manipulating the activities of Cln2 or Cdc20, they found that [Clb2] can be either in high or 
low, depending on initial conditions. In terms of bifurcation theory, they provided evidence 
for an S shaped steady state curve bounded by saddle-node bifurcations, with transitions 
driven by the activity of Cln2 or Cdc20. Indeed, we found that Chen's model displays such 
bifurcations when [Cln2] and [Cdc20] are considered as bifurcation parameters. 

In accord with the experimental protocol of Cross et al. Q], we consider [Cln2] and 
[Cdc20] as parameters, and therefore discard Eqn. (Al) and Eqn. (A4-A5) from Chen's 
model. We performed bifurcation analysis for the remaining six ODE's. In Figure 3, we 
show a combination of two bifurcation diagrams. In the left bifurcation diagram we set 
[Cln2]=0 and vary [Cdc20], whereas in the right bifurcation diagram we set [Cdc20]=0 and 
vary [Cln2]. As mass is the same in both cases (M = 1), the two stable steady states in 
Figure 3 represent Gl phase (low Clb2-dependent kinase activity) and S/G2/M phase (high 
Clb2-dependent kinase activity). Figure 3 shows that increasing [Cln2] drives the transition 
from Gl into S/G2/M, while activation of [Cdc20] drives the transition from S/G2/M back 
to Gl. 

Using AUTO's facility for computing two-parameter bifurcation diagrams, we extended 
the saddle-node bifurcations in Figure 3 into the parameter planes spanned by ([Cln2],M), 
([Cdc20],M), and ([Cln2],[Cdc20]). In Figure 4a, there are multiple steady states inside the 
cusp-shaped region bounded by the dashed lines, as expected 1^. In Figure 4b, there are 
two different bistable domains, bounded by dashed and dotted hues, respectively. Where 
the domains overlap, we found that the control system has five steady states. We found 
that two different modules independently lead to the bistable domains in Figure 4b. The 
dashed line curve is due to the Hctl module of the wiring diagram in Figure 1, whereas the 
dotted line curve is due to the Sicl module. Finally, Figure 4c shows the bistable region on 
the ([Cln2],[Cdc20]) plane. 
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V. EFFECTS OF TRANSCRIPTIONAL FACTORS MCMl AND SBF IN A RE- 
DUCED MODEL WITH FOUR ODE'S 



Because Eqn. (A1-A20) take into account many known details of cell cycle control, the 
model is very complex. It is difficult to understand from Eqn. (A1-A20) what are the 
nonlinearities leading to specific features of the bifurcation diagram shown in Figure 2. To 
overcome this difficulty, we simplify Chen's model, by defining a core module that retains 
the main dynamical features of the full set of equations. The reduced model can be useful 
in understanding the roles of nonlinear feedbacks in the control system. 

In Figure 5 we propose a simplified wiring diagram for the budding yeast cell cycle. We 
discarded from the original wiring diagram the Sicl and the Clb5 modules and Cdc20's 
activation, retaining only four ODE's. 



-[Cln2] = M(fc;„2 + <n2[SBF]) - fc,,„2[Cln2], 

|[Clb2] = M{K^,, + ^.jMcml]) - {k',^,, + {k'l,, - k',^,,mctl] + fc:,%[Cdc20])[Clb2], 

l.^^^^. ^ (fcUi + fe:.a[Cdc20])(l-[Hctl]) _ y,a[Hctl] 
dt^ ^ J,,a + 1 - [Hctl] J,,ii + [Hctl] ' 

4cdc20] = (A;;2o + <2o[Clb2]) - k',^,,[Cdc20], 



where [SBF] is given by Eqn. (A13-A14) with [Clb5]=0. [Mcml] is given by Eqn. (All), 
and Vi^ti is given by Eqn. (A15). With the elimination of the dynamics for Cdc20 activation, 
we define a new parameter in Eqn. (4) k'^ 20 = k 2o^+v' 20+fcd 20 ' ^ote that the results in 
this section do not change if k'^20 — ^d,2o- 

Although Figure 5 is much simpler than Figure 1, Eqn. (1-4) are still quite complex. 
The most uncertainties arise from the two transcription factors(SBF and Mcml), which are 
described by nonlinear Goldbeter-Koshland functions 3)^^- Their role is to switch solu- 
tions from one branch to another. As the effects of the transcription factors can be studied 
experimentally, we explore their roles via two-parameter bifurcation diagrams. First, by 
using mass as the primary bifurcation parameter, we computed a one-parameter bifurcation 
diagram similar to Figure 2. Then, we continued the codimension-one bifurcations into two 



parameter domains, using the intensity coefficients of the transcription factors, k'!. ,^2 
k'!.^J2, as the secondary bifurcation parameters. 

Figure 6a shows a two-parameter bifurcation diagram of Eqn. (1-4) on the (M, 
plane. Saddle-node bifurcations appear in Figure 6a only if ^s„2 > 0- ^ bistable domain 
is inside the dashed lines. It widens at smaller mass and larger fc^'„2- The dot-dashed line 
showing Hopf bifurcation points continues inside the bistable domain. At M > 0.8, the Hopf 
bifurcation line is accompanied by a cyclic fold curve. These curves eventually coalesce at a 
larger mass value. Inside the bistable domain the cyclic fold coalesces with a locus of saddle 
loops (solid line in Figure 6a). 

Figure 6b shows a two-parameter bifurcation diagram on the (M, A;s,b2) plane. A crucial 
difference between Figure 6a and Figure 6b is the existence of a bistable domain at fc"^2 — 0- 
If k'g i,2 < 0.04, the effect of [Mcml] regulation is negligible. But if fc"j,2 > 0.5, [Mcml] can 
destroy bistability. In Figure 6b a Hopf bifurcation line originates from a Bogdanov-Takens 
bifurcation. This line is accompanied by a line of saddle loops. The saddle loops change 
stability where the line of cyclic folds coalesces with the line of saddle loops. 

We found that Eqn. (A1-A20) display two-parameter bifurcation diagrams similar to 
Figure 6a-b. Notice from Figure 6b that the domain of bistability is quite independent of 
the activity of Mcml, but the existence of the primary Hopf bifurcation in the model is 
sensitively dependent on the activity of Mcml. 



VI. A SNIC BIFURCATION IN A REDUCED MODEL WITH THREE ODE'S 



The eukaryotic cell cycle engine is a highly conserved molecular machine. It is expected 
that mathematical models of cell cycle controls in different organisms exhibit qualitatively 
similar dynamics as revealed by similar bifurcation diagrams. But there can be also pecu- 
liarities in these models, subject to particular parameter selections. As we mentioned in 
Section 3, the bifurcation diagram in Figure 2 does not involve a SNIC bifurcation, as seen 
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in bifurcation diagrams of mathematical models for frog eggs and fission yeast 
Although this difference is rather subtle and does not contradict any features of cell cycle 
physiology, we point out that Chen's model Eqn. (A1-A20) can display a SNIC bifurcation 
for appropriate choice of parameter values (not shown). In this section, we examine SNIC 
bifurcation in a three variable model. 
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To further simplify the model, we neglect Cln2 from the wiring diagram in Figure 5. As 
a result, we have a model with three time-dependent variables. 



-[Clb2] = M(K^,, + kl,,[M<^ml]) - (k',^,, + (kl,, - fc^,,2)[Hctl] + A;:;;,,[Cdc20])[Clb2], 

. ^ «a + Ci[Cdc20])(l-[Hctl]) _ \^.,i[Hctl] 
dt^ ^ J„,n + 1 - [Hctl] J,,a + [Hctl]' 

icdc20] = (fc^^o + <2o[Clb2]) - fc,,2o[Cdc20]. 



In Eqn.(5-7), [Mcml] is given by Eqn. (All), and Vi^ti is given by Eqn. (A15). We 

k' M 

assume [Cln2] — and [Clb5] = — in Eqn. (A15). We also changed the values of some 

d,b5 

parameters in Table I, as fc^ (,5 = 0.06, k"^^^ = 1.5, k"2o = 0-07, Ja,mcm = Ji,mcm = 0.01. 

Let [Clb2]°, [Hctl]° and [Cdc20]° denote a steady state solution of Eqn. (5-7). Clearly, 
^C^c20f = K2o+kUic^^^f) _ Substituting this functional relation between [Cdc20] and [Clb2] 
into Eqn. (5-6), we can think of ([Clb2],[Hctl]) as a two- variable system, susceptible to phase 
plane analysis. The nuUclines of the two variable systems are plotted in Figure 7. From 
the intersections of these nuUclines, we find steady state solutions, [Clb2]° and [Hctl]°, 
and consequently [Cdc20]°. Depending on M, the number of intersections varies, but the 
maximum number of steady states is three. 

We study stability of the steady states numerically. In Figure 8 we plot a bifurcation 
diagram of Eqn. (5-7), with M as the principal bifurcation parameter. At a given M, there 
can be one stable steady state and two unstable steady states, or a single steady state which 
can be either stable or unstable. The stable steady states in Figure 8 can coexist with stable 
limit cycle oscillations. There are two interesting features in this bifurcation diagram: (i) 
a SNIC bifurcation which arises at M 2.9 where a saddle- node coalescence is replaced 
by limit cycle oscillations, and (ii) byrhithmicity, i.e., coexistence of two stable limit cycle 
oscillations, for 2.8 < M < 4.9. 

Figure 9 shows a two-parameter bifurcation diagram for Eqn. (5-7). Despite the reduction 
to just three ODE's, this diagram is quite complex. Multiple steady states are found inside 
the sohd hues. There are three Bogdanov-Takens bifurcations in Figure 9, from which 
originate three independent loci of Hopf bifurcations, shown by lines in violet. Two cyclic 
folds associated with the Hopf bifurcations are shown by lines in cyan. Two saddle loops. 
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shown by green lines, originate at BTi and BT2, cross the region of bistabihty, and attach 
to the right saddle node line at two saddle-node-loop bifurcation points. Between these two 
saddle-node-loops, we find a SNIC bifurcation, red line in Figure 9. 



VII. BIRHYTHMICITY AND EFFECTS OF EXTRINSIC FLUCTUATIONS 

Figure 9 shows that the distance between the inner cyclic fold(CF2) and the SNIC bi- 
furcation varies as changes. In other words, depending on k'!. fj2^ birhythmicity may 
occur either close to, or far from the START transition (when the stable Gl state gives 
way to large amplitude stable oscillations). If it happens far away from START, it will not 
interfere with cell cycle progression. However, if it occurs close to START, as in Figure 8, 
an interesting question arises. To which stable oscillation (the large amplitude or the small 
amplitude limit cycles) will the trajectory of motion (see red lines in Figure 2) connect? 
We found that the trajectory of motion always follows the large amplitude slow oscillations 
in the three-variable model. We have shown (in a separate publication) that switching be- 
tween small and large amplitude oscillations is possible when the model takes into account 

r~i 

diffusion terms |24]. Here, we demonstrate the effects of noise on the trajectory of motion. 

A complex process, such as cell cycle control, is naturally subject to fluctuations from 
different sources. For instance, stochastic effects due to size and nuclear volume differences 
at cell division have been studied for fission yeast Since we know very little about 

the origin of fluctuations in the cell cycle engine, the simplest way to incorporate random 
processes into Eqn. (5-7) is to assume that certain extrinsic fluctuations randomly perturb 
the cell cycle engine. Mathematically, we replace Eqn. (5-7) by Langevin-type equations 



with multiplicative noise 



26 



27| 



^ [Clb2] = F[cib2] + pD,[C\h2]i{t) , (8) 

^[Hctl] = F[Hcti] + v/2I^2[Hctl]e(t), (9) 

^[Cdc20] = F[cdc20] + 72D3[Cdc20]e(t), (10) 

— = ,M. (11) 
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where, F[cib2], -^[Hcti], -P'[cdc20] are the right hand sides of Eqn. (5-7) and ^(t) is Gaussian 
white noise with zero mean and unit variance, 



<e(t) = o>, <^(t)at')>=sit-t'). (12) 



We assume that mass increase is not affected by random fluctuations 



We simulated Eqn. (8-11) usine standard numerical techniques for stochastic differential 

nn 

equations [2^ |2^. In Figure 10 we overplot two different simulations. The dashed lines 
show time evolutions of M, [Clb2], [Cdc20] and [Hctl] when birhythmicity occurs far from 
START. In this case, noise does not interfere with mitosis, and cell mass divides each time 
[Clb2] drops below 0.1. The solid lines show the case when birhythmicity occurs close to 
START, as in Figure 8. In this case noise can switch the control system from slow, large 
amplitude oscillations to fast, small amplitude oscillations. As a result, [Clb2] does not go 
below 0.1 and the cell cannot divide. Consequently, mass M grows and the system goes 
to the stable steady state (see filled diamonds at M > 3.9 in Figure 8). Therefore, in the 
presence of noise, birhythmicity may lead to mitotic arrest. 



VIII. DISCUSSION 



In this work, we carried out bifurcation analysis of a model of the budding yeast cell 
cycle, based on earlier work by Chen et al. which successfully accounts for many observed 
features of proliferating yeast cells. Our results show that, despite a pecuharity in topology 
of the bifurcation diagram, the budding yeast cell cycle model displays the same basic 
features previously associated with frog egg and fission yeast models; namely, saddle-node 
bifurcations associated with stable and unstable oscillations. 

We explored bistability and hysteresis in this model by numerical bifurcation analysis. 
Some of our bifurcation diagrams can be useful for designing new experiments. For instance, 
our two parameter bifurcation analysis (Figure 4b) suggests that the [Hctl] and [Sicl] mod- 
ules may lead independently to bistable states, and there can be regions in parameter space 
with three stable steady states, when these two modules operate cooperatively. 

We found that a reduced model with four time- dependent variables retains the main char- 
acteristics of the bifurcation diagram of Chen's model. This reduction allows us to explore 
the dominant roles of SBF and Mcml transcription factors in budding yeast checkpoint 
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controls. Our two-parameter bifurcation diagrams (Figure 6) also can be useful in designing 
experiments for cell cycle controls by transcription factors. 

The budding yeast cell cycle model of Chen et al. is parameter rich. Although the 
parameter set presented in Table I leads to a satisfactory fit of the model to many exper- 
imental observations, the choice of parameter values should be further constrained by new 
biochemical data about the protein-protein interactions and further improved by automatic 



parameter estimation techniques 



On the other hand, different sets of parameters. 



leading to different bifurcation scenarios, are interesting from a theoretical standpoint. We 
have proposed a set of parameters for a reduced, three-variable model leading to a SNIC 
bifurcation. 

An interesting feature accompanying the appearance of a SNIC bifurcation in the reduced 
model is birhythmicity. Birhythmicity has been found in a chemical system js^, but for 



15, lla, 133 ■ We have shown that in the 



biological systems, it is known theoretically only 
presence of extrinsic fluctuations, birhythmicity can lead to mitotic arrest. The fact that 
noise can switch a biochemical system from one stable solution to another is well known 
(e.g. Ref. QQ), 

but switching from one stable oscillations to another is a less studied 
research area. A more systematic study of switching between stable limit cycles is a problem 
for the future. 
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Figure Captions 
Fig. 1. Wiring diagram of a budding yeast cell cycle model j^. 

Fig. 2. A one-parameter bifurcation diagram of Eqn. (A1-A20) for parameter values 
in Table I. Solid lines indicate stable steady states. Dashed lines indicate unstable steady 
states. Solid circles denote the maximum and minimum values of [Clb2] y on stable limit cycle 
oscillations, open circles denote the same for unstable oscillations. The red line shows the 
trajectory of motion when Eqn. (A1-A20) are supplemented by the mass growth equation 
M = /iM. The cell divides {M = f ■ M) when [Clb2]T drops below 0.005. 

Fig. 3. Bistability and hysteresis driven by [Cln2] and [Cdc20]. On the left plane 
[Cln2] = 0, on the right plane [Cdc20] = 0. Mass is fixed at 1. Filled diamonds show stable 
steady states, dashed lines show unstable steady states. Dotted lines and arrows indicate 
the START and FINISH transitions of the hysteresis loop. (START refers to the Gl — > S 
transition, FINISH refers to the M — > Gl transition.) 

Fig. 4a. Two-parameter bifurcation diagram on the ([Cln2],M) plane. Multiple steady 
states are found inside the cusp-shaped curve. 

Fig. 4b. Two-parameter bifurcation diagram on the ([Cdc20],M) plane. There are two 
independent pairs of saddle-node bifurcation curves in this figure (dashed curves and dotted 
curves). Depending on the overlaps of the regions bounded by these curves, the number of 
steady states varies from I to V. 

Fig. 4c. Two parameter bifurcation diagram on the ([Cln2], [Cdc20]) plane. Bistable 
steady states are found in between the dashed curves. 

Fig. 5. Wiring diagram of a reduced model with four- variables. 

Fig. 6a. Two-parameter bifurcation diagram of Eqn. (1-4) on the (A;"„2) ^) plane. Loci 
of saddle node bifurcations (SN) are shown by dashed lines. Other lines trace saddle loop 
(solid), cyclic fold (dotted), and Hopf (dot-dash) bifurcation points. 

Fig. 6b. Two-parameter bifurcation diagram of Eqn. (1-4) on the {k'!.^2i^) plane. 
Bistable steady states are found inside the SN curve (dashed line). Other lines indicate loci 
of saddle loops (solid), cyclic folds (dotted), and Hopf bifurcations (dot-dash). The locus of 
Hopf bifurcations originates from a Bogdanov-Takens point shown by the filled circle. 

Fig. 7. Stationary solutions of Eqn. (5-6) can be computed from the intersections of the 
[Hctl] nuUcline (solid line) and the [Clb2] nullcline (dashed line). In this plot mass is fixed 
at M = 2. Notice that [Clb2] 0.15 is the region where [Mcml] changes abruptly from 
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to 1. 

Fig. 8. Bifurcation diagram of Eqn. (5-7). Filled diamonds show stable steady states, 
dashed hnes show unstable steady states. Stable limit cycle oscillations are shown by filled 
circles, unstable limit cycle oscillations are shown by open circles. 

Fig. 9. Two-parameter bifurcation diagram of Eqn. (5-7) on the (M, /c^'^g) plane. Bista- 
bility is found inside the SN curve (solid line). Three different Hopf bifurcations (violet 
lines) originate from three Bogdanov-Takens bifurcation points shown by filled circles at 
BTi, BT2 and BT3. Cychc folds are shown by lines in cyan, saddle loops by lines in green, 
and the red solid line shows SNIC bifurcations. SL3, which runs next to HB3, is not shown 
on the diagram. CF2 runs from a degenerate Hopf bifurcation on HBi to a degenerate Hopf 
bifurcation on HB2. CFi runs from a degenerate saddle loops on SLi to a degenerate saddle 
loop on SL2 (not shown), crossing over HBi on the way. Where CFi and CF2 run very close 
together, only CF2 is plotted on the figure. 

Fig. 10. Stochastic simulations of Eqn. (8-11). Dashed hnes show a case when birhyth- 
micity occurs far from the SNIC bifurcation. In this case, noise does not interfere with cell 
cycle progression. Solid lines show the case when birhythmicity occurs close to the SNIC 
bifurcation, as in Figure 8. In the presence of noise, the latter case leads eventually to 
mitotic arrest. Parameters are: Di — D2 — D3 — 3.75 • 10~^. Solid lines for kg ,j2 = 0.05, 
dashed hnes for k"fj2 — 0.06. 
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APPENDIX A: MODEL FOR BUDDING YEAST CELL CYCLE 



|[Cln2] = M(A;;„2 + KM^BF]) - ^,,„2[Cln2], (Al) 

|[Clb2]T = M(A;;,2 + <.2[Mcml]) - (fc^,,^ + (^d,.2 - ^d,.2)[Hctl] + k'^^,,[Cdc20])[C\h2]T, (A2) 
d . . _ K,n + <a[Cdc20])(l - [Hctl]) l^,,i[Hctl] 

c^t ^ ^ Ja,i + 1 - [Hctl] + [Hctl] ' ^^"^ 

^[Cdc20]T = (/c^^o + <2o[Clb2]) - fc,,2o[Cdc20]T, (A4) 

^[Cdc20] = A;„,2o([Cdc20]r - [Cdc20]) - (V^,,2o + A;d,2o)[Cdc20], (A5) 

|[Clb5]T = M(A;;,5 + fc^'^^jMBFj) - {k',^,, + A;:;,,5[Cdc20])[Clb5]T, (A6) 

4sicl]T = A:;,, + <,JSwi5] - (A:,i,ei + , ^'l.^ )[Sicl]T, (A7) 

at Jd2,cl + [olCijT 

^[Clb5|Sicl] = fc„,,b5[Clb5][Sicl]- 

-{kdi,b, + ^^,65 + fcd,65[Cdc20] + + , ^'t 11 )[Clb5|Sicl], (A8) 

^[ab2|Sicl] = A;„,,b2[Clb2][Sicl] - 

-{kdi,m + (A:^,,2 + {Km - fcd,62)[Hctl] + fcX62[Cdc20]) + k,,,,, + )[Clb2|Sicl], (A9) 

[Cln3]* + eei,fe2[Bck2] + [Cln2] + eei,65[Clb5] + 6^,62 [Clb2]), (AlO) 

[Mcml] — G(fc(j^j^cm[Clb2] , ki^^fncrm Ja,mcmi Ji,mcm)i (■'^-^■^) 

[Swi5] = G(A;„,,^,[Cdc20], k^^, + A:,%jClb2], J„,,^,, J,,,^,), (A12) 

[SBF] = [MBF] = GiVa,sbf, k^.f + <,,^[Clb2], J„,,,/, J,,,,y), (A13) 

K,s6/ = ^a,s6/([Cln2] + esbf,n2, ([Cln3]* + [Bck2]) + e,;./,65[Clb5]), (A14) 

Vi,ti = + <a([Cln3]* + e,,n,n2[Cln2] + Q,n,65[Clb5] + e,,a,62[Clb2]), (A15) 

[Clb2]T = [Clb2] + [Clb2|Sicl], (A16) 

[Clb5]T = [Clb5] + [Clb5|Sicl], (A17) 

[Sicl]T = [Sicl] + [Clb2|Sicl] + [Clb5|Sicl], (A18) 

[Bck2] = M[Bck2]°, (A19) 

[Cln3]* = [Cln3]^„^ — — . (A20) 
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TABLE I: Kinetik constants for the budding yeast model 



Rate constants (min ^] 



k' = 0.006 



K,n2 — 



K,b5 = 0-02 
k'' = 0.05 



b2 



0.002 



k'd,b2 = 0-010 
K,ci = 0.020 
kas,b2 = 50 
K,20 = 0.005 

ka,2Q = 1 



h" 

h" 



0.05 

2 



<ci = 0.1 

kas,b5 — 50 

fc^o = 0.06 



k' 



20 



a.tl 



ka,sbf 



k' 



i,sbf 



0.04 



0.5 



h" 

ka,mc 



0.1 

= 2 



^:^,65 = 0.1 

kd,n2 =0.1 



^dV = 0-05 
kdi,ci = 0.01 
fcdi,62 = 0.05 
kd,20 = 0.08 
10 



Ksbf = 6 



'^1,20 

ka,su 
k'- 







1 

0.3 



Characteristic concentrations (dimensionless) 
[Cln3]max = 0.02 [Bck2]0 = 0.0027 

Ja,sbf ~ Ji,sbf ~ 0.01 Ja,mcm ~ Ji,mcm ~ 1 

Ja,swi ~ Ji,swi ~ 0.1 Ja,tl ~ Ji,tl ~ 0.05 

Kinase efficiencies (dimensionless) 



Ccl.nS = 20 
^i,tl,n2 = 1 
^sbf,n3 = 75 



^cl,k2 = 2 
^i,tl,b2 = 1 
es6/,b5 = 0.5 



Other parameters 

/ = 0.433 J„3 = 6 



eci,62 = 0.067 
ej,ii,fe5 = 0.5 
Vi,2o = 0.05 

£'n3 = l 



^^5 = 0.25 



kd2,ci = 0.3 
kdi,b5 = 0.05 



ki,mc 
k'! 



0.64 
= 0.15 
= 0.2 



J, 



d2,cl 



0.05 



ecl,65 



/X = 0.005776 



GoUdbeter Koshland function: 



G{a, b, c, d) = '^'"^ = (A21) 

h-a-\-bc-\-ad-\-^{h-a-\-bc-\- ad^ - Aab{h - a) 
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